HDL-CR-89-001-1 


March  1939 


U-Azi 


Solution  of  Two-Dimensional,  Two-Region  Electromagnetic 
Ground  Response 


ty  ira  Konlberg 


DTIC 

ELECTE 

APR  1  2198S 


Prepared  by 

Kohiberg  Associates,  Inc. 
F.O.  Box  23077 
Alexandria,  VA  22304 

Under  contract 
DAA 103-38-0-0001 


U.S.  Army  Laboratory  Command 
Harry  Diamond  Laboratories 
Adeiphi,  MD  20783-1197 


Approved  for  public  release;  distribution  unlimited. 


fl  T  ■  O 


*  1 


REPORT  DOCUMENTATION  PAGE 


Form  Approved 
OMP  No.  0704-0188 


1*.  REPORT  SECURITY  CLASSIFICATION 

UNCLASSIFIED 


2a.  SECURITY  CLASSIFICATION  AUTHORITY 


2b.  DECLASSIFICATION /DOWNGRADING  SCHEDULE 


1b  RESTRICTIVE  MARKINGS 


3  DISTRIBUTION /AVAILABILITY  OF  REPORT 


Approved  for  public  release;  distribution  unlimited. 


4  PERFORMING  ORGANIZATION  REPORT  NUMBER(S) 


6b.  OFFICE  SYMBOL 
(if  tppliable ) 


TCN  87-116 


6a.  NAME  OF  PERFORMING  ORGANIZATION 
Kohlberg  Associates 


6c  ADDRESS  (City,  Stitt,  end  ZIP  Code) 

P.O.  Box  23077 
Alexandria,  VA  22304 


8a.  NAME  OF  FUNDING  /SPONSORING 
ORGANIZATION 

Harry  Diamond  Laboratories 


8c.  ADDRESS  (City,  Stttt,  tnd  ZIP  Code) 

2800  Powder  Mill  Road 
Adelphi,  Maryland  20783-1197 


1 1  TITLE  (Include  Security  dissipation) 

Solution  of  Two-Dimensional,  Two-Region  Electromagnetic 


8b  OFFICE  SYMBOL 
(if  tppliable) 

SLCHD-NW-EP 


5.  MONITORING  ORGANIZATION  REPORT  NUMBER(S) 

HDL-CR-89-001-1 


7a  NAME  OF  MONITORING  ORGANIZATION 

U.S.  Army  Research  Office 


7b.  ADDRESS  (City,  Sate,  end  ZIP  Code) 

PO  Box  12211 

Research  Triangle  Park,  NC  27709-2211 


9.  PROCUREMENT  INSTRUMENT  IDENTIFICATION  NUMBER 


10  SOURCE  OF  FUNDING  NUMBERS 


PROGRAM  PROJECT  TASK 

ELEMENT  NO  NO1L162120  NO 

6.21. 20A  AH25 


Ground  Response 


13a.  TYPE  OF  REPORT  13b  TIME  COVERED  _  _  _  _ 

Final  from  Mar  87  to  Aug  87  I  March  1982 


14.  DATE  OF  REPORT  (Yeir,  Month,  Dey) 

March  1989 


15.  PAGE  COUNT 
39 


12.  PERSONAL  AUTHOR(S) 
Ira  Kohlberg 


16.  SUPPLEMENTARY  NOTATION 

This  task  was  performed  under  a  Scientific  Services  agreement  issued  by  Battelle,  Research  Triangle  Park  Office,  200 
Park  Drive,  P.O.  Box  12297,  Research  Triangle  Park,  NC  27709.  HDL  Project  No.  XE77E3,  AMS  Code  612120.H25001 1 


17. _ COSATI  CODES _  18.  SUBJECT  TERMS  ( Continue  on  reverse  if  necessity  tnd  identify  by  block  number) 

field  group _ sub-group  EMP,  electromagnetic  pulse^ Maxwell's  equations,' boundary  condition',  Green's 

09 _ 03 _ function,  multilayer  ground,  earth  model;  HABEMP,  high-altitude  burst  EMP 

20  14  -/  'Lri  ^ 


19.  ABSTRACT  (Continue  on  reverse  if  necessity  tnd  identify  by  block  number) 


~~~A  solution  for  the  two-dimensional,  two-region  electromagnetic  ground  response  has  been  developed  which  relates 
the  surface  components  of  the  electric  field  to  the  surface  components  of  the  magnetic  field.  This  has  been  accomplished 
by  deriving  a  universal  functional  form  for  a  dimensionless  Green's  function.  The  Green's  function  provides  increasingly 
more  accurate  approximations  to  the  response  for  each  successive  reflection  from  the  second  layer.  This  result  would 
appear  to  provide  simplification  and  reduced  computer  running  time  in  the  numerical  modelling  of  the  HABEMP  when  the 
ground  response  is  coupled  to  finite-difference  methods  for  solving  the  atmospheric  part  of  the  problem. 


20  DISTRIBUTION /AVAILABILITY  OF  ABSTRACT 
IS  UNCLASSIFIED/UNLIMITED  □  SAME  AS  RPT 


22*.  NAME  OF  RESPONSIBLE  INDIVIDUAL 
William  T.  Wyatt 


DO  Form  1473,  JUN  86 


121  ABSTRACT  SECURITY  CLASSIFICATION 

UNCLASSIFIED 


22b  TELEPHONE  (Include  ArtiCode)  22c  OFFICE  SYMBOL 
(703)490-2303  SLCHD-NW-EP 


Previous  editions  ire  obsolete 


SECURITY  CLASSIFICATION  OF  THIS  PAGE 

UNCLASSIFIED 


CONTENTS 

1.  INTRODUCTION  1 

2.  MATHEMATICAL  NOTATION  2 

3.  EQUATIONS  AT  SURFACE  OF  EARTH  6 

4.  MATHEMATICAL  STRUCTURE  OF  Gn{x,t)  11 

5.  CALCULATION  OF  Gn(x,t)  17 

5.1  Calculation  of  20 

5.2  Calculation  of  Git2(x,t)  20 

5.3  Calculation  of  G\t3(x,t)  22 

6.  CONCLUSION  28 

ACKNOWLEDGEMENT  29 

APPENDIX  A  —  CALCULATION  OF  Q(k,  t)  33 

APPENDIX  B  —  APPLICATION  OF  FALTIJNG  THEOREM  35 

DISTRIBUTION  37 

FIGURES 

1.  (n  +  l)-region  model  3 

2.  Geometric  considerations  3 

3.  Ho (td,4>d)  as  a  function  of  <f>d  with  td  as  parameter  for  r  =  0.5  32 

4.  n0(trf,^)  as  a  function  of  t &  with  4>d  as  parameter  for  r  =  0.5  32 


TABLES 

1.  Values  of  IIq for  r  =  0.3,  0.5,  0.7,  2.0,  5.0 


1.  INTRODUCTION 


In  a  previous  study1  we  showed  that  for  an  n-layer  earth  model  it  is  possible  to  express 
the  three  components  of  the  electric  field  and  the  vertical  component  of  the  magnetic  field 
on  the  surface  of  the  earth  as  a  space-time  integration  of  the  two  horizontal  components 
of  the  magnetic  field.  In  particular,  it  was  shown  that  if  is  a  point  on  the  surface  ( x ,  y 
plane) ,  and  if  Y,-  is  a  member  of  the  set 

Yi{ra,t)={Hz,Ex,Ey,Ez}  ,  (1.1) 

then  every  member  of  the  set  Y,  on  the  surface  of  a  finitely  conducting  earth  can  be  related 
to  the  horizontal  components  of  the  surface  magnetic  field  through  the  equation 

Yx{ra,t)  =  f  f  Gix(fa  -  r\,t  -  t')Hx{rl,t')  dx'  dy'  dt' 

JO  J rjj 

+  /  /  Giv[f*s  -  r'9,t  -  t')Hy{r'a,t')  dx'  dy'  dt'  .  (1.2) 

Jo  J r* 

The  functions  G,x  and  G,v  are  Green’s  functions,  which  are  determined  from  the  solution 
of  the  n- region  ground  model  shown  in  figure  1. 

The  result  given  by  equation  (1.2)  can  possibly  provide  considerable  simplification  in 
the  numerical  modelling  of  the  high-altitude  burst  electromagnetic  pulse  (HABEMP)  when 
the  ground  response  is  coupled  to  finite-difference  methods  for  solving  the  atmospheric  part 
of  the  problem.  When  this  approach  can  be  used  it  obviates  the  necessity  of  developing  a 
numerical  representation  of  the  ground,  which  then  reduces  the  number  of  variables  in  the 
problem  and  hence  the  computer  running  time  (and  cost) .  On  the  other  hand,  the  reduction 
of  machine  variables  must  be  weighed  against  the  speed  of  the  numerical  computation  for 
»  the  integral  boundary  conditions  arising  from  the  Green’s  function  formalism.  This  ques- 

U)l.  Kohlberg,  Surface  Integral  Representation  of  Three  Dimensional  Electromagnet¬ 
ic  Ground  Response  for  Multi-Layered  Earth  With  Frequency  Dependent  Electrical  Param¬ 
eters,  April  1986,  Final  Report  for  Harry  Diamond  Laboratories,  Contract  No.  DAAG29- 
81-D-010,  Delivery  Order  2064,  Battelle  Columbus  Laboratories  (Prime  Contractor). 
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tion  is  as  yet  unresolved;  in  any  event  it  will  depend  on  the  number  of  ground  layers  being 
considered,  the  values  of  the  electrical  parameters,  and  the  time  range  of  interest. 

Equation  (1.2)  is  derived  by  solving  the  »region  three-dimensional  problem  in  the 
Fourier  (space)  and  Laplace  (time)  domains  and  ultimately  performing  inverse  Fourier/La¬ 
place  transforms.  (The  reader  is  referred  to  reference  1  for  details.)  Although  the  formalism 
and  mathematical  procedure  are  general,  the  analytical  Fourier/Laplace  transform  inver¬ 
sion  may  not  always  be  possible.  For  the  purposes  of  assessing  the  feasibility  of  using 
the  theory  of  reference  1,  four  special-case  solutions  were  developed.  These  included  the 
following  models: 

(a)  one  layer,  one  dimension, 

(b)  one  layer,  two  dimensions, 

(c)  one  layer,  three  dimensions,  and 

(d)  two  layers,  one  dimension. 

In  this  study  we  extend  the  number  of  solutions  to  include  the  two-layer  two-di¬ 
mensional  case.  The  geometry  for  this  case  is  shown  in  figure  2.  Following  an  earlier 
analysis,1  we  consider  the  case  where  the  displacement  current  in  the  ground  is  less  than 
the  conduction  current.  We  assume  that  spatial  variations  in  the  y-direction  are  neglected. 
Our  principal  concern  is  to  establish  the  space-time  relationships  between  Ex,Ey  and  the 
driving  functions,  Hx,Hy. 

2.  MATHEMATICAL  NOTATION 

As  will  become  evident  in  sections  3  and  4  of  this  report,  several  mathematical  oper¬ 
ations  are  performed  on  the  four  surface  values  of  the  fields  used  in  this  analysis:  namely, 
Ex(x,t),  Ey(x,t),  Hx(x,t),  and  Hy(x,t).  The  purpose  of  this  section  is  to  define  these  op¬ 
erations  and  develop  a  shorthand  notation  for  dealing  with  them.  If  F(x,t )  represents  any 
one  of  the  aforementioned  four  functions,  the  Fourier  and  Laplace  transforms  are  defined 
by  the  following  operations: 


i  r°° 

F{k,t)  =  MF(x,t )  =  — -  /  e~xkxF{x,  t)  dx 

27T  J_oq  ■ 


(2.1) 


2 


z-n  +  l  ►  OO  on  +  i,en+unn+1 

Figure  1.  (n  +  l)-region  model. 
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Figure  2.  Geometric  considerations. 


F{x,s) 


LF(X,t)  =  r 
J  o 


e  stF(x,t)  dt 


M  =  Fourier  transform  operator,  and 
L  =  Laplace  transform  operator. 
The  double-transformed  function,  F(k,s ),  is  given  by 


(2.2) 


F(k,s )  =  M_LF(x,t)  =  LMF(x,t) 

1  /*  oo  /•  oo 

=  —  j  e-5 *  J  e~'kxF(x,t)  dx  dt 


(2.3) 


As  indicated  in  equation  (2.3),  the  order  in  which  the  transforms  are  taken  is  imma¬ 
terial. 

The  inverse  Fourier  and  Laplace  operators  are  defined  by  operations 


/oo 

eikxF(k,t)  dk  =  F(x,t)  , 

-OO 


-j  rC  +  lOO 

L~1F{x,s)  =  —  /  estF(x,s)  ds  =  F(x,t) 

Z'Kl  J c  —  ioo 


It  also  follows  that 


(2.4) 


(2.5) 


M  F(fc,s)  =  F{x,s)  , 


(2.6) 


L~lF{k,s)  =  F{k,t)  .  (2.7) 

Using  the  Faltung  theorem  we  can  write,  for  any  functions  f\{k,s)  and  f2{k,s), 


1  . 

M  1  {fi{k,s)f2(k,s))  =  —  J  fi(x  -  x',s)f2{x',s)  dx' 

=  —  /  /2(x  -  x',s)/i(x',s)  dx'  .  (2.8) 

27r  ./-oo 
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The  shorthand  notation  for  the  operation  of  equation  (2.8)  is 


M  1  (/i(M)/2(M))  =  h{x,s)  <2>  h{x,s)  ,  (2.9) 

where  the  “®”  stands  for  the  space  integration  of  equation  (2.8).  It  also  readily  follows 
that 

M-1  (/i(M)/*(M))  =  ®  /2(x,t)  •  (2.10) 

The  convolution  theorem  for  the  Laplace  transform  yields  the  following  result: 


L~ 1  (/i(M/2(M)  =  f  /i(M  ~  t')f2(k,t')  dt' 

J  o 

=  f  A(M  -  t')fi  {k,t')  dt' 

Jo 


(2.11) 


The  shorthand  notation  for  the  convolution  integration  in  equation  (2.11)  is 


LT1  (/i(M/2(M)  =  /i(M)  *  h{k,t)  ,  (2.12) 

where  denotes  the  convolution  operation.  We  also  have 

L~l  (/i(^,s)/2(2:,s))  =  fi[x,t)  *  f2{x,t)  .  (2.13) 

The  mathematical  formalism  and  notation  developed  through  equation  (2.13)  provides 
a  compact  way  of  identifying  the  Green’s  function  which  relates  the  surface  values  of  the 
electric  field  to  the  magnetic  field  components.  For  example,  as  shown  in  section  3,  we 
have  the  relationship 

Ey(k,s)  =  G(k,s)Hx(k,s)  .  (2.14) 

Using  equations  (2.8)  and  (2.9)  we  can  write 

M-1  L~lEy{k,s)  =  Ey(k,t)  =  Ey(x,t)  .  (2.15) 
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If  we  identify  G(k,s)  as  fi(k,s)  and  Hx(k,s )  as  /2(Ar,s),  we  deduce  the  relationship 


M~x  L~1Ey{k,s)  =  M_1  L~x  ( G{k,s)Hx{k,s )) 
Ey(x,t )  =  G(x,t)  0  *Hx(x,t ) 

*oo  r  t 


poo  pi 

Ey(x,t)  =  /  /  G(x  —  x' ,t  —  t')Hx(x' ,t')  dx'  dt'  (2.16) 

J  —  oo  J  0 


where 


G(x,t)  =  L~1G(k,s)  =  L~x  M~xG{k,s) 

—  Green’s  function  .  (2.17) 

Employing  the  techniques  of  reference  1  we  can  derive  an  expression  for  G(k,s).  The 
thrust  of  this  investigation  is  to  develop  techniques  for  determining  G[x,t)  from  the  inverse 
Fourier/Laplace  transforms. 


3.  EQUATIONS  AT  SURFACE  OF  EARTH 

It  is  shown  in  equation  (3.21)  of  reference  1  that  for  the  nth  region  of  an  n-region 
ground  model,  the  equations  for  EXiTl(z,kx,ky,s)  and  Ey<n{z,  kx,  ky,  s)  in  terms  of  the 
magnetic  field  components,  Hx,n{z,kx,ky,s),  Hyin{z,kx,ky,s),  are  given  by 


EXyn  — 


In 


mHy<n  ~f~  2  {kykxHz,n  + 

In 


(3-1) 


E  —  — 

&y,n  ~ 

u  n 


where  on  —  conductivity  and 


Hx,n  o  {k XH x<n  "I"  kxkyHy}Tl) 
In 


Tn  —  i  An 


(3.2) 


(3.3) 


An  —  \j Sfln(7n  ~f“  ^ y  ’ 


(3.4) 


with  +  being  used  for  the  upward  wave,  and  —  being  used  for  the  downward  wave. 
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We  immediately  note  from  equations  (3.1)  and  (3.2)  that  if  either  ky  or  kx  is  equal 
to  zero  (this  corresponds  to  neglecting  spatial  variations  in  the  y  and  x  directions,  re¬ 
spectively)  the  equations  are  decoupled  into  two  independent  pairs.  For  example,  setting 
ky  —  0  gives 


ET>  _  _ tt 

^1,11  —  ^t/,n  > 

G  n 


£y,n 


V  'll' 


(3.5) 


(3.6) 


with  '/n  now  being  given  by 


In  =  ±\/sflnC Tn  +  kl  . 


(3.7) 


When  equations  (3.5)  to  (3.7)  are  used  in  the  solution  of  the  two-region  problem  of 
figure  2,  it  can  be  shown  that  the  surface  field  equations  are  given  by 


Ex  =  rx 


Hy  +  2  E 


n=  1 


Ey  =  —W\ 


Hx  +  2^W 


where 


Ai 

r  l  =  —  , 
o\ 

Sjl 

Ml  =  T"  » 


n  =  l 


A2 

r  2  =  — 
02 


(3.8) 


(3.9) 


w2  = 


SfJ, 

Ao 


Ai  •=  \fSfio i  +  k2  ,  A2  =  v sfxo2  +  k2  , 
r2  -  r  1 


= 

*«,  = 


r2  +  rx 
w2  —  wx 


exp(— 2Aj/)  , 

exp(— 2A  i /)  . 


(3.10) 


W2  +  W\ 

For  brevity  we  have  replaced  kx  by  k,  it  being  understood  that  we  only  considering 
spatial  variations  in  the  x  direction. 
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Equations  (3.8)  and  (3.9)  reduce  to  the  one-layer,  two-dimensional  case  considered  in 
reference  1  when  <72  =  cr\-  In  this  situation,  'I',,  and  are  both  equal  to  zero,  and  we 
obtain  the  solutions 

Ex,0  =  r1Hy  ,  (3.11) 

Ey,  o  =  -u»i  Hx  .  (3.12) 

The  space-time  behavior  of  Eyio{x,  t)  is  given  by1 


Ey,0  A /  L  Eyfi 


H  f*  1 
27T  Jo  J t  -  t'exp 


/  t\  2 

( X  —  X  )  \lO  i 

4(f-i')  j 


dx'  dt'  .  (3.13) 


For  the  one-layer,  two-dimensional  case,  equation  (3.13)  represents  a  means  of  es¬ 
tablishing  the  boundary  condition  on  top  of  the  earth,  ^his  equation  can  be  written  in 
finite-difference  form  and  thus  can  be  used  with  the  HABEMP  equations  above  the  earth’s 
surface  to  provide  a  self-consistent  representation  of  the  overall  physical  model.  The  key 
question  is  “Under  what  conditions  is  the  method  more  efficient  than  representing  the 
ground  in  a  finite-difference  approximation?”  The  answer  to  this  will  be  forthcoming  in 
the  near  future. 

For  the  two-layer,  two-dimensional  case  being  considered  here,  the  deduction  of  Ey(x,  t ) 
in  terms  of  Hx(x,t)  is  more  complicated  than  that  of  equation  (3.13)  because  of  the  com¬ 
plexity  of  s).  A  similar  statement  can  be  made  for  the  Ex,Hy  pair;  however,  for 

brevity  this  pair  is  not  being  considered  since  the  analysis  is  similar  to  the  Ey,Hx  case. 

It  is  possible  to  express  the  solution  of  equation  (3.9)  in  several  different  ways,  leading 
in  turn  to  different  algorithms  for  completing  the  calculation  of  Ey(x,t).  For  example, 
performing  the  L-1  M-1  operation  directly  on  equation  (3.9)  leads  to  the  result 


OO 

Ey  =  Eyt 0  +  2  Eyyn  , 

n=  1 


(3.14) 
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where 


Ey<n  =  Gn(x,t )  ®  *Hx(x,t)  ,  n  >  1  ,  (3.15) 

Gn(x,t)  =  L~lM'1  (-Kwt)  .  (3.16) 

An  alternative  and  mathematically  equivalent  formalism  is  based  on  using  the  knowl¬ 
edge  of  Ey< o  deduced  from  equation  (3.13).  We  can  write  equation  (3.9)  in  the  form 

OO 

Ey  =  Eyt0  +  2j2KEyio  ,  (3.17) 

n—  1 

which  then  leads  to  the  relationship 

OO 

Ey  ~  Ey>0  +  2  ^2  Ey>n  •  (3.18) 

n  =  l 

Eytn  is  expressed  as 

Ey,n  =  Gn(x,t)®*Ey,o(x,t)  ,  (3.19) 

with  Gn(x,t )  now  given  by 

Gn(x,t)  =  L~lM~1  (**)  •  (3.20) 

The  analysis  of  this  investigation  is  concerned  with  the  determination  of  Gn(x,t)  from 
equation  (3.19). 

There  is,  however,  one  additional  method  that  can  be  used,  which  is  related  to  equation 
(3.20).  We  include  this  for  completeness.  Consider,  for  example,  the  sequence  of  functions 
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1 


Eytl  —  tymEy'Q 

Ey,2  =  **^,0  =  *wEy}0 


Ey,n  ^  wEy<n—i 

Using  equation  (3.21)  we  can  write 

OO 

Ey  —  Ey> o  +  2  Ey<n 

n=  1 

E y >n  is  now  given  by 


Ey,n  —  G  i  <S>  *Eytn-i 

with 


(3.21) 


(3.22) 


(3.23) 


C?i(x,<)  =  '&w(x,t)  =  L  1 M  1'&w(k,s)  .  (3.24) 

It  is  also  observed  by  comparing  equations  (3.20)  and  (3.24)  that  Gn(x,t )  is  the 
(n  —  1)3<  space-time  convolution  of  ^’lu(x,t).  That  is, 


G2  =  L-1M-1{*l)=*w®**w  , 
G$  =  <8>  *&w  ®  * &w  , 


n  —  1  convolutions 

a  —  vif  ,  — ^ —  -s  vi/ 

un  —  *  ^  io 


(3.25) 
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In  summary,  it  is  clear  that  the  determination  of  'F w(x,t ),  as  given  by  equation  (3.24), 
is  the  basic  building  block  of  the  calculation.  This  is  the  focus  of  the  effort  of  the  next 
section. 


4.  MATHEMATICAL  STRUCTURE  OF  Gn(x,t) 
From  equation  (3.10)  we  can  write 

=  {9{k,s)v[k,s))n  =  gnvn  , 

where 


9n{k,s) 


W2  —  U>  i 
W2  +  W  i 


(4.1) 


(4.2) 


vn(k,  s)  —  exp(-2nAj/)  .  (4.3) 

Examination  of  equation  (3.17)  shows  that  the  electric  field  at  the  surface  can  be 
considered  as  a  sum  of  terms  involving  multiple  round-trip  reflections  from  the  second 
layer.  This  can  be  seen  for  example  by  first  examining  the  inverse  Laplace  transform  of 
vn{k,s). 

We  have 


vn(k,t)  =  L  1  vn(k,  <s) 

=  L~1exp(—2nl\/sfiai  +  k 2) 

=  L_:1exp(  —  2nlyJixo\ y/s  +  Qi)  ,  (4.4) 


where 


k2 

a  i  =  - 

/ZCTl 


Using  the  formula 


(4.5) 
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L~‘f(s+ at)  = 


(4.6) 


we  obtain 


vn(k,t)  =  e  ait6n{t)  ,  (4.7) 

where 


and 


MO 


2s/n  t 3/2 


exp 


~n2TH\ 

4*  y 


(4.8) 


T/?  =  L2  no j  =  two-way  diffusion  time  to  second  layer  , 

L  =  21  —  round-trip  distance 

As  observed  from  equation  (4.8),  for 

~  >  1  (4.9) 

the  damping  will  become  severe  and  perhaps  only  one  term  in  the  expansion  will  be 
necessary.  Moreover,  6n  is  in  general  a  rapidly  decreasing  function  of  n,  as  can  be  seen  by 
examining  its  maximum  value,  On,ma x-  This  is  determined  by  solving  the  equation 


d.0n 

dt 


=  0 


for  the  time  at  which  the  maximum  occurs;  this  time  is  given  by 


t 


m  — 


n2TR 

6 


and  the  corresponding  value  of  0n  is 


63/2  1  1  _a 

2y/n  TR  n2 


(4.10) 


(4.11) 


(4.12) 
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The  foregoing  equation  supports  the  conjecture  that  higher-order  terms  provide  diminish¬ 
ing  contributions  to  the  overall  solution. 


The  main  difficulty  in  determining 


Gn(x,t)  =  L~l  M'1*”  = 


(4.13) 


is  attributed  to  performing  the  inverse  Laplace  transform  of  gn{k,  <s).  Using  equation  (3.10) 
in  equation  (4.2)  yields 


gn(k,s)  =  ^ 


Al  ~  A2 
X\  +  A2 


(4.14) 


where 


Ai-Aj 


i  9  . 

=  xmw  (41S) 

and  A i ,  A2  are  given  by  equation  (3.10).  Let  us  first  consider  some  of  the  mathematical 
properties  of  g*.  We  write  gi  in  the  form 


Ai  —  A2  Aj  —  Aj 

Ai  -f-  A2  Ai  —  Aj 


(Aj  —  A2)2 

M -A! 


(4.16) 


Using  equation  (3.10)  we  have 


A?  -  \l  =  ~  . 


(4.17) 


Substituting  equations  (4.16)  and  (4.17)  into  equation  (4.14),  and  subsequently  using  the 
binomial  expansion  for  (Aj  —  A2)2ri  we  obtain 


/  1  \  n  2n 

=  E  ^r~<-A,r  . 


(4.18) 


where  f>2n,m  is  the  binomial  coefficient 


^2  n,m  — 


(2  n)! 


(4.19) 


(Zn  —  m)\m\ 

Let  us  now  consider  the  structure  of  the  terms  in  equation  (4.18).  For  the  even  terms, 
characterized  by 


m  =  2m'  ,  0  <  m’  <  n 


(4.20) 


we  have 


A2n-2m'(_A2)2m'  =  ^  +  +  *2)”*'  .  (4.21) 

For  the  odd  numbered  terms,  characterized  by 

m  =  2m  +1  ,  0  <  m  <  n  —  1  ,  (4.22) 

we  can  write 


A?n-m"(-A2)m"  =  -{spun  +k2)n  m"  1  {sfio2  +  k2)m"  sT  (4.23) 

where 


-  _  y/sn<7 1  +  fc2y/s/ig2  +  fc2 
s 

By  expressing  (s/zax  +  fc2)n  m  ,  [spo2  +  fc2)m  ,  (spgj  +  /c2)m  ,  and  (s/zo2  +  k2)m 
as  a  binomial  expansion,  and  substituting  the  results  into  equation  (4.18),  it  can  be  shown 
by  combining  the  summations  that  gn  is  of  the  form 


(4.24) 


9n  —  9n  1  d"  9n2  > 


where 


9n  1  — 


1 


Sn{n{<7 1  -  02)) 


<r= 0 


rk2p 


(4.25) 


9n  2 


Sn{p{o  1  -  cr2))’ 


^2  Brsrk2q  sT 


(4.26) 


In  the  foregoing  expressions,  Ar  and  Br  are  constants  which  depend  on  ax,a2  and  the 
index  r,  and  the  integers  p  and  q  are  positive  linear  functions  of  the  summation  index,  r. 
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It  is  not  necessary  to  go  into  the  tedious  details  to  establish  the  implications  of  equations 
(4.25)  and  (4.26)  regarding  the  computation  of  Gn(x,t)  from  equation  (4.13). 

Using  equations  (4.1)  and  (4.3)  in  equations  (4.25)  and  (4.26)  we  have 

\T/n  —  \T/rl  _L  ■5/n 
*10  *to,l  '  *10,2 

where 


1  /  fr2p  \ 

n>i  =  ^Ar  Vn(k's})  ’  (4-27) 
1  (  k 2q  _\ 

■  <4-28) 

From  equation  (4.13)  we  deduce 


where 


If  we  now  let 


Gn(x,t )  =  Gn\  (x,t)  +  Gn2  (x,  t )  , 


1 

G  nl  =  /  \  \  n  ^  ^  ArYr  (x,  t)  , 

(M<u  ~  ^2))  v 


n-  1 


^  ~  M°i  -  °2))n  ^  BrZ,(x,t)  ’ 


/  fc2p  \ 

Yr{x,t )  =  M"1  L"1  (  —  vn(/c,s)J  , 

(  k2q  _\ 

Zr(x,t)  =  M  lL  1  vn{k,s)sT J 


and 


Vn(z,0  =  M  1  L  1Vn(fc,-s) 


(4.29) 

(4.30) 

(4.31) 

(4.32) 

(4.33) 

(4.34) 
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$n(z,0=M  1L  1  [vn(k,  s)sr(fc,s))  ,  (4.35) 

it  then  follows  from  equations  (2.6)  and  (2.7)  that 


Yr(x,t)  =  j  •••  J 

f  dt  i  dt%  •  dtfi — f  | 

f  i 

(  3x0  Mx'h)  ’ 

(  d2  \  q 

(4.36) 

1 

Zr(x,t)  =  J  •••  J 

'  dt  1  dt 2  •  •  •  dtn-r  ( 

k  3xi)  ' 

(4.37) 

In  summary,  we  have  shown  that  if  one  can  determine  vn(x,t)  and  <£n(x,  £),  it  is  pos¬ 
sible  to  determine  the  Green’s  function  Gn(x,t)  by  space  derivates  on  vn  and  followed 
by  repeated  time  integrations.  The  utility  of  this  approach  depends  on  the  simplicity  and 
speed  of  performing  the  time  integration. 

In  equation  (4.7)  we  showed 

vn(k,t)  =  L~1vn{k,s)  =  exp  t)  0n(t)  .  (4.38) 

V  J 

We  now  have 

vn(x,t)  =  M~1  L~1vn(k,s)  =  0n(t)R(x,t)  ,  (4.39) 

where 


R(x,t)  =  J  elkx  exp  ^ 

(  x2noi 

—  — -  exp 


..  ,  ,  •  (4.40) 

Vi  \  U  J 

It  is  observed  that  vn(x,t)  can  be  expressed  in  closed  form,  which  facilitates  the  compu¬ 
tation  of  equation  (4.36). 

Now  let  us  consider  the  deduction  of  <3>n(x,£).  We  have 


$»(*><)=  M  1  L  1  (svn{k,s)T(k,s)) 
=  §-t{M-1L~1Pn(k,s ))  , 


(4.41) 
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where 


Pn{k,s)  =  vn(k,s)T(k,s)  .  (4.42) 

By  examining  equations  (4.7),  (4.8),  and  (4.24),  which  are  the  constitutents  of  equa¬ 
tion  (4.42),  we  notice  that  the  only  difference  between  Pn  and  P\(k,s)  is  the  replacement 
of  Tr  by  n2TR.  Thus,  if  we  can  determine 


Pi(x,t)  =  M-1  L~1Pi  (k,  s) 


(4.43) 


we  can  determine  Pn(x,t)  by  substituting  u2Tr  for  Tr  in  the  resultant  expression. 

For  the  purposes  of  this  investigation  we  are  limiting  the  calculation  to  the  evaluation 
of 


Gi{x,t)  =  M  1  L  1  (gi(k,s)vi(k,s))  .  (4.44) 

This  provides  the  contribution  to  the  surface  value  of  the  electric  field  from  the  first  round- 
trip  reflection.  In  addition,  section  5  shows  that  the  calculation  of  Gi(x,t )  includes  as  one 
of  its  components  the  computation  of  JPi(x,t). 

5.  CALCULATION  OF  G t(x,t) 

Starting  from  equation  (4.44)  we  have 


Gi(x,t)  =  M-1  L~1(g1(k,s)v1(k,s)) 

=  MT1  [  gi{k,t')vi(k,t  —  t')  dt' 

Jo 

/OO  ft 

etkx  /  gi(k,t')vi(k,t  —  t')  dt ’  dk 

-oo  Jo 


Using  equation  (4.16)  we  can  write 


9 1  = 


sii{oi  -  a 2) 


(Aj  -f  Aj  —  2Aj A2)  , 


(5.1) 


(5.2) 
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•S, 


(5.3) 


Si  = 


V  x  "  V  '  L  V  1  '■*■4 


where 


sfi{o i  -  ct2) 

k 2 


<*i  = 


a2  = 


M<7i 

fc2 

MCT2 


We  also  have 


«  /,  °\  +  °2  r/ *  2 k2  rr/  .  2v/ai<72  ,  - 

Si(M)  =  - - —6(t)  +  — - —H{t)  -  —  L~l Q 


<7i  —  <72 


~  02) 


0\  —  o2 


where 


Q  =  V5  +  ai  Vs  +  q2  _  1  p 

■S  llsfO\02 


S(t)  =  delta  function  , 
H(t)  =  step  function  . 

From  equation  (5.5)  we  can  write 


(5.4) 


(5.5) 


(5.6) 


9i[k,t)  —  <7n  +  pi2  +  g13  , 


where 


9n  -  - - —  S(t)  , 


9 12  = 


°  1  ~  CT2 
2k2 

n(ai  -  a2) 


m  . 


9 13  =  <?(M)  , 

0\  —  02 


Q(k 


,t)  =  L~'  ( 


\/s"  +  Oil  ■y/s  +  a2 


As  shown  in  appendix  A,  Q(k,t)  is  given  by 


(5.7) 


(5.8) 


Q(M)  =  e-a<  [/?/i(/?f) +<*/<,(/?*)] 

+  (a2  -  P2)  f  e~au I0((3u)  du  +  6{t) 
Jo 


(5.9)  | 

I 
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where 


Using  equations  (4.7)  and  (4.8)  we  can  write 

«,(M)=ex  p(-^)(5^)«p(-^)  (5.11) 

If  we  insert  equations  (5.7)  to  (5.11)  into  equation  (5.1),  and  note  that  we  can  combine 
the  8{t)  of  equation  (5.9)  with  that  of  equation  (5.8),  we  can  write 


G'i(I>  t)  =  G\ti  +  Gi'2  +  Gi,3  ,  (5.12) 

where 


oo  rt 


G  ij  =  /  etkx  /  4>j(k,t')v\(k,  t  —  t')  dt'  dk  ,  j  = 
J  —  oo  J  0 

1,2,3  ,  (5.13) 

(5.14) 

<y  Jl2 

(5.15) 

4>z {k,t)  -  '  ^  \Li(k,t)  +  L2{k,t)  +  L3{k,t) 

0\  —  <72  L  -1 

(5.16) 

Ll(k,t)=Pe-atI1((3t)  , 

(5.17) 

L2{k,t)  =ae~atI0{Pt)  , 

(5.18) 

L3(k,t)  =  (a2  -  /?2)  [  e~auI0(Pu)  du  . 

Jo 

(5.19) 
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It  should  be  noted  that  although  o\  —  a 2  appears  in  the  denominator  of  equations 
(5.15)  and  (5.16),  the  sum  of  <^2  +  4>z  must  equal  zero  in  the  limit  of  o\  =  (note  that 
<f>  1  is  already  zero  in  this  case).  This  follows  from  the  observation  that  there  can  be  no 
reflections  in  this  situation.  This  is  easily  seen  by  using  approximations  to  Lj,L2,and  L3 
in  the  limit  of  (3  — *  0. 

In  the  next  subsections  we  perform  the  calculations  of  Gij. 

5.1  Calculation  of  Ghl(x,t) 

From  equations  (5.13)  and  (5.14)  we  have 


[  <j)1(k,t)vi(k,t  —  t')  dt'  =  [  6(t')vi  (k,  t  —  t')  dt' 

Jo  \f°  1  +  Jo 

__  y/57  ~  y/02 


k2t 

. —  . :exp  1 - 

y/Ol  +  v/^2  V 


where 


Performing  the  integration  over  k  space  yields 


fit) 


,u\  V^R  { -Tr\ 

m  =  2^  eXP \~S~ )  • 


(5.20) 


(5.21) 


G  -  yfil-  Vgj 
1,1  yj°l  +  y/^2 


mi 


eikxexp 


k( 


HO  1 


dk 


=  Vor-yg?  x  / 

yfo 7  +  y/oi  \  At  )  y/i 

_  y/dT  ~  \/^2  1  Tr  1 

~  v^T+V^  2  7^  Iexp 


Tr  /.  x2\ 


5.2  Calculation  of  G1>2 

From  equations  (5.13)  and  (5.15)  we  have 


/oo  <■  t 

e'kl  / 

-00  ./  0 


2k2 


Hi°i  ~  02) 
Replacing  k2  by  -d2  jdx2  we  obtain 


exp 


k2f 

Ho  1 


f{t')  dt'  dk 


(5.22) 


(5.23) 
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G\,2  —  - 


f)  2  ft  roo 

a?Lmit'Ljk^ 


dk 


v(<n  -  02)  a r2  •/-«,  '  “'*■  W,  j  -"  (5-24) 

The  integral  in  equation  (5.24)  is  the  same  as  that  in  equation  (5.22),  so  that  we  can  write 


k2t' 


Gl  2  — 


a2 


1 ' "  /x(<7l  —  CT2)  dx 2 

Making  the  substitution 


pi  JTr_  1 

Jo  2  (t')2  L 


exp 


1  Tr(\  +  fi) 


(5.25) 


u  4  ^  V  +  X2  J  *' 


(5.26) 


permits  the  integration  to  be  performed  through  the  relationship 


du'  — - 7T  dt 


We  obtain 


/  2 


=  -7T 


(5.27) 


a2 


/iL  o\  —  (7 2  ax2 


i/(x,0  , 


(5.28) 


where 


U  (x,  t)  =  — - — —  exp  (  -  Tr ^  +  L 2  ^ 

1  +  jj 


4  *  ; 

Using  equation  (3.19),  the  contribution  from  the  first  reflection  will  be  given  by 


(5.29) 


Eyti  —  Gi  (x,t)  ®  *Eyto(x,t)  , 
and  in  particular  the  contribution  from  G  1,2  is 


(5.30) 


EyA'2~*L  oi-o2l_00l0(dx'*Ui<X,,t')')Ey'o{x  X>,t  t^dx'dt'  ■  (5<31) 

Based  on  application  of  the  Faltung  theorem  discussed  in  appendix  B,  the  foregoing 
integration  can  be  converted  to  the  form 
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E 


V.1,2  - 


nL  o\ 


*-■  •  (5-32) 


It  may  turn  out  that  for  computational  purposes  the  form  given  by  equation  (5.32)  is  easier 
to  evaluate. 


5.3  Calculation  of  G1i3 

Inserting  equation  (5.16)  into  equation  (5.13)  and  rearranging  gives  the  following 
expression  for  G  i)3: 


Gl)3(I)0  =  _W  dt'  f(t-t,)(Nl(x,t')  +  N2(x,t')  +  N3(x,t'))  ,  (5.33) 

0\  —  02  Jo 


where 


Ni  =  [  fielkxe  at  li(0t')exp  (—  — — (t  —  t')^  dk  ,  (5.34) 

J- oo  V  J 

N2=[0°aeikxe-at'l0{/3t')exp(-  —  (t-t,))dk  ,  (5.35) 

J- oo  \  J 


N$  =  [  (a2  —  P2)elkx  [  e  auI0((3u)  du  exp  — — (f  —  t')^  dfc 

»  —  oo  J  0  V  / 

equations  (5.34)  to  (5.36)  can  be  simplified  using  the  formulas 


r°(2) = £  r 

^  7o 


e*cos*d0  , 


w=l-[ 

n  Jo 

and  writing  a  and  /?  in  the  form 


eZCO30  cos  0d&  , 


(5.36) 


(5.37) 


(5.38) 


a= - A:2  , 


(5.39) 
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(3=  - k 2  , 

Hob 


where 


(5.40) 


1  _  1  /  1  1  \  _  2(<71ct2) 

o  l  ’ 

cra  2  \CTi  cr2  /  +  <?1 

1  _  1  (  1  1  \  _  _  2(<t1ct2) 

/>  I  )  > 

°b  2  \£7i  ct2/  cr2  —  <T X 

Using  equations  (5.37)  to  (5.42)  in  equations  (5.34)  and  (5.35)  yields 


1  1/1 


(5.41) 


(5.42) 


where 


,r  a2 
N‘  =  “ai2  n‘  ' 
a2 

^  =  “ai2  "2  ’ 


x  exp  (  (- — )k2t'  cos  6  )  exp 

V  f**b  / 


.  ^°l  / 


(5.43) 


rii  =  — - —  [  cosO  dd  [  e‘fezexp  ( —  )k2t,'\ 

V-W  Jo  J-o o  \  / 

X  exp  f (— —  )/cV cos 0^  exp  f  — — — (i  —  t1)]  dk  ,  (5.44) 

\  flOb  )  \  i  J 

n2  =  _J_  f  do  f°°  eikxexp  (~(—)k2t') 

Jo  J-oo  \  ) 


(5.45) 


The  expression  for  N3  is  simplified  by  first  performing  the  integration  over  u  in  equa¬ 
tion  (5.36).  We  have 


/  e~auI0{(3u)  du  =  -  f  dO  l  e~auepucos6  du 
Jo  n  Jo  Jo 

l  f*  i  _,e—  (q— P  cos  e)t' 

—  —  I  d6 - 

7 r  J0  a  —  (3  cos  0 


(5.46) 
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Inserting  equation  (5.46)  into  equation(5.36)  yields 


where 


JV3  = 


a2 

dx*n3 


n3  -  -i—  [  dO  S{0 )  [  eikxexp  ( — 

[10aTT  JO  J  —  oo  \ 


A:2(t  —  t') 


[1(7 1 


X  (l  -  exp  -  — )*V))  , 

V  \  /  ) 


S(0)  = 


1  -  (^): 

v  CTfc  / 


1  “  (ft)cos  0 


Using  the  general  formula 


/: 


■* = ^ exp  \ « 


Akx-Ak*  -  '  ~X 

e  e 


=  ^(A)  , 


with 


A,  =  —  fl-^cosfl>\«,+  —  (<~*#)  - 

MCTa  V  °b  J  [K7 1 


A2 


f  - 

/i<7l 


we  deduce 

n3  =  — C"F(A2) - —  f  d9  S{0)F{ Ax)  , 

[10aK  [l(7an  Jo 

where 


Noting  that 


(7  a 

02  —  O  i 

Ob 

cr2  +  0\ 

(5.47) 


(5.48) 


(5.49) 


(5.50) 


(5.51) 


(5.52) 


(5.53) 


(5.54) 


(5.55) 
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gives 


and 


dd 

1  -  ( cra/ob )  cos  6 


7 r 

sjl  -  (oa/afc)2 


(5.56) 


C"  =  ?r\/l  -  ( oalob )2  (5.57) 

By  comparing  equation  (5.51)  with  the  exponential  terms  of  equations  (5.44)  and 
(5.45)  we  obtain 

1 

rii  —  -  /  cos  $  ^0  F(A.i)  ,  (5.58) 

y0 


1  /■*■ 

n2  =  -  /  d0  F(A2) 

noan  J  o 

From  equations  (5.53),  (5.58),  and  (5.59)  we  have 


(5.59) 


n  =  rei  +  n2  +  ^3 


(5.60) 


v/l  -  {oa/ob)2  1  [* 

n  =  — - F(A2)  + -  /  fF(0)F(Ax)  d0 

pcra7r  y0 


where 


(5.61) 


w(o)  =  —  cos  « + 1  -  s(e)  , 

Ob 

iyW=  KM)^in20  , 

1  —  ( Oa/Ob )  COS  0 


(5.62) 


Combining  equations  (5.43),  (5.47),  (5.60),  (5.61)  and  inserting  the  latter  into  (5.33)  gives 


where 


G 


1,3  - 


dx2 


[. ABh  +  ACh) 


(5.63) 
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A  _  2 y/OlOl 


Cf  1  —  02 


B  = 


C  = 


\J  1  -  KM-)2  x/i  -  p2 


(gg/gfe) 

llOa'n 


P-O* 
2 


H<7a 


lie  an 


h  =  f'dt'  f(t-t')F( A2) 
Jo 

h=  f  dt'  f(t-t')(F{ Ai)) 
Jo 


(F(Al))=  f  T^p li-eF{Al) 

Jo  l  —  p  cos  0 


p=  —  = 


pc  OS  0 
CT2  —  O’  1 


(5.64) 


tffc  0-2  +  0 1 

It  is  easy  to  see  that  Ji  reduces  to  the  same  integral  as  that  of  equation  (5.25);  that 


is, 


h  =  £ 

where  U(x,t)  is  given  by  equation  (5.29). 

Using  equations  (5.41)  and  (5.42)  we  deduce 


(5.65) 


d2(ABIi)  4  d2  ^  ,  . 

~ap~  =  smw['(i',)  = 


(5.66) 


We  observe  that  the  foregoing  contribution  provides  an  exact  cancellation  to  Glt 2. 

Thus,  the  significant  contribution  to  Gii3  comes  from  the  term  ACI2.  Unfortunately, 
this  does  not  appear  to  be  easily  integrated,  and  numerical  techniques  must  be  used. 
For  simplicity  in  presentation  we  shall  render  the  results  in  dimensionless  form.  We  first 
compute  the  constant  AC;  we  have 


AC  = 


o\  —  o2 


K 


Pyj0\02  <J\  +  <72 


(5.67) 


where 
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K  = 


/Oi  (7 1  —  <72 
a2  al  +  °2 


Let  us  now  introduce  the  variable 


V 

V=1 

in  equation  (5.64d).  We  can  write  I2  in  terms  of  y  in  the  following  form: 


I2  — 


where 


n”  =  ^/„’(r^exp(W^)) 

<t>d' 


n{<t>d,y)  -  J  y~~ 


lo  =  1  + 

T2 

4>d  =  ^- 
4  td 

t 

U~TR 

x 

L 
02 


0  (1-1/)* 

K  sin2  0  d.0  1 

p  COS  6 
1  -  r 


exp  - 


n{<t>d,y)  dy  , 


lo 


2  r 


t/(l  +  COS0) 


Xd~  L 


r  — 


P  = 


o  1 

r  -  1 
r  +  1 


Using  the  dimensionless  variables  introduced  in  equation  (5.71)  we  write,  for  the 
function, 


Gi  =  G  1,1  +  ^2  (ACI2) 


(5.68) 


(5.69) 


(5.70) 


(5.71) 

Green’s 


(5.72) 
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Gx 


_1  1-y/r  1  1  (  l  +  xg\ 

2  i  +  trl  qexp  v  ; 

1  _1_  1  -  r  1  d 2 

2  v/r  1  —  r  TrL  dx\  0 


(5.73) 


The  first  term  in  equation  (5.73)  is  that  of  G1(1  (cf.  equation  (5.22c)),  expressed  in  terms 
of  Tr,L,  and  the  dimensionless  variables  r,  t d,  and  Xd. 

The  contribution  to  the  electric  field  from  Gi  is  given  from  equation  (5.30),  where 
the  required  integrations  can  be  performed  in  either  Xd/td  or  x/t  space.  It  is  also  evident 
from  equations  (5. 31)  and  (5.32)  that  the  term  involving  d2W.0jdx\  can  be  converted  to 
an  integration  of  II0  combined  with  the  second  partial  derivative  of  £y,0. 

Since  IT0  cannot  be  expressed  in  closed  form,  it  must  be  provided  in  tabular  form, 
either  as  a  function  of  td,  Xd  or  td ,  <f>d-  This  consideration  should  be  reserved  when  the  actual 
implementation  of  this  result  is  incorporated  in  the  HABEMP  finite-difference  model. 

Tabular  values  of  Tl0{td,<{>d)  for  selected  values  of  td  and  (f>d  are  rendered  in  table  1 
(pp  30  and  31).  Plots  of  II0  as  a  function  of  <f>d  with  td  as  a  parameter  for  r  =  0.5  are 
shown  in  figure  3,  while  plots  of  II0  as  a  function  of  td  with  ^  as  a  parameter  for  r  =  0.5 
are  shown  in  figure  4. 


6.  CONCLUSION 

A  solution  for  the  two  region,  two-dimensional  electromagnetic  ground  response  has 
been  developed  which  relates  the  surface  components  of  the  electric  field  to  the  surface 
components  of  the  magnetic  field.  This  has  been  accomplished  by  deriving  a  universal 
functional  form  for  a  dimensionless  Green’s  function.  The  Green’s  function  provides  in¬ 
creasingly  more  accurate  approximations  to  the  response  for  each  successive  reflection  from 
the  second  layer.  This  result  would  appear  to  provide  simplification  and  reduced  computer 
running  time  in  the  numerical  modelling  of  the  HABEMP  when  the  ground  response  is 
coupled  to  finite-difference  methods  for  solving  the  atmospheric  part  of  the  problem. 
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Table  1.  Values  of  Tlo{td,4>d)  for  r  =  0.3,  0.5,  0.7,  2.0,  and  5.0 


r  =  0.3 


4>d 

0.00 

0.25 

0.50 

0.75 

1.00 

1.50 

2.00 

3.00 

4.00 

5.00 

0.10 

0.4189 

0.3450 

0.2847 

0.2354 

0.1949 

0.1345 

0.0935 

0.0463 

0.0237 

0.0126 

0.30 

1.7276 

1.4584 

1.2341 

1.0469 

0.8902 

0.6485 

0.4771 

0.2657 

0.1533 

0.0915 

0.50 

2.0805 

1.7740 

1.5166 

1.2999 

1.1171 

0.8312 

0.6247 

0.3630 

0.2184 

0.1355 

0.75 

2.1579 

1.8532 

1.5957 

1.3776 

1.1925 

0.9003 

0.6865 

0.4104 

0.2537 

0.1615 

td 

1.00 

2.1258 

1.8340 

1.5865 

1.3760 

1.1965 

0.9117 

0.7015 

0.4268 

0.2682 

0.1734 

1.50 

1.9976 

1.7334 

1.5081 

1.3156 

1.1506 

0.8867 

0.6899 

0.4288 

0.2750 

0.1810 

2.00 

1.8634 

1.6227 

1.4169 

1.2459 

1.0951 

0.9649 

0.7545 

0.5954 

0.3801 

0.2497 

3.00 

1.6254 

1.4212 

1.2459 

1.0951 

0.9649 

0.7545 

0.5954 

0.3801 

0.2497 

0.1681 

4.00 

1.4368 

1.2592 

1.1065 

0.9747 

0.8607 

0.6760 

0.5356 

0.3446 

0.2280 

0.1544 

5.00 

1.2952 

1.1370 

1.0006 

0.8828 

0.7808 

0.6151 

0.4888 

0.3161 

0.2102 

0.1430 

r  =  0.5 


4>d 

0.00 

0.25 

0.50 

0.75 

1.00 

1.50 

2.00 

3.00 

4.00 

5.00 

0.10 

0.4091 

0.3340 

0.2730 

0.2235 

0.1832 

0.1236 

0.0839 

0.0393 

0.0189 

0.0093 

0.30 

1.7094 

1.4274 

1.1942 

1.0010 

0.8407 

0.5965 

0.4266 

0.2233 

0.1205 

0.0669 

0.50 

2.0696 

1.7447 

1.4738 

1.2476 

1.0584 

0.7664 

0.3054 

0.1718 

0.0993 

0.75 

2.1547 

1.8288 

1.5556 

1.3260 

1.1328 

0.8319 

0.6160 

0.3459 

0.2000 

0.1187 

td 

1.00 

2.1278 

1.8140 

1.5498 

1.3270 

1.1387 

0.8438 

0.6305 

0.3603 

0.2119 

0.1278 

1.50 

2.0055 

1.7193 

1.4772 

1.2720 

1.0977 

0.8226 

0.6215 

0.3630 

0.2179 

0.1340 

2.00 

1.8744 

1.6125 

1.3902 

1.2012 

1.0401 

0.7848 

0.5969 

0.3531 

0.2146 

0.1335 

3.00 

1.6385 

1.4151 

1.2249 

1.0625 

0.9236 

0.7022 

0.5381 

0.3229 

0.1989 

0.1253 

4.00 

1.4502 

1.2553 

1.0890 

0.9467 

0.8248 

0.6298 

0.4846 

0.2931 

0.1818 

0.1153 

5.00 

1.3084 

1.1343 

0.9856 

0.8582 

0.7488 

0.5735 

0.4426 

0.2692 

0.1679 

0.1069 

r  =  0.7 


4>d 

0.00 

0.25 

0.50 

0.75 

1.00 

1.50 

2.00 

3.00 

4.00 

5.00 

0.10 

0.4063 

0.3301 

0.2685 

0.2187 

0.1783 

0.1189 

0.0797 

0.0363 

0.0169 

0.0080 

0.30 

1.7093 

1.4189 

1.1797 

0.9824 

0.8195 

0.5732 

0.4036 

0.2042 

0.1061 

0.0565 

0.50 

2.0754 

1.7385 

1.4590 

1.2266 

1.0331 

0.7368 

0.5292 

0.2789 

0.1509 

0.0837 

0.75 

2.1651 

1.8257 

1.5423 

1.3054 

1.1069 

0.8004 

0.5829 

0.3158 

0.1757 

0.1001 

td 

1.00 

2.1409 

1.8130 

1.5382 

1.3075 

1.1136 

0.8123 

0.5968 

0.3290 

0.1862 

0.1078 

1.50 

2.0213 

1.7210 

1.4682 

1.2549 

1.0746 

0.7925 

0.5887 

0.3316 

0.1915 

0.1131 

2.00 

1.8911 

1.6155 

1.3828 

1.1859 

1.0190 

0.7565 

0.5657 

0.3227 

0.1887 

0.1127 

3.00 

1.6550 

1.4193 

1.2196 

1.0500 

0.9056 

0.6774 

0.5103 

0.2953 

0.1750 

0.1059 

4.00 

1.4658 

1.2598 

1.0849 

0.9360 

0.8091 

0.6078 

0.4597 

0.2682 

0.1601 

0.0975 

5.00 

1.3231 

1.1389 

0.9823 

0.8488 

0.7348 

0.5536 

0.4200 

0.2463 

0.1478 

0.0905 
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Table  1.  Values  of  n0 (tdi<f>d)  for  r  =  0.3,  0.5,  0.7,  2.0,  and  5.0 

(continued) 

r  =  2.0 


4>d 

0.00 

0.25 

0.50 

0.75 

1.00 

1.50 

2.00 

3.00 

4.00 

5.00 

0.10 

0.4258 

0.3428 

0.2762 

0.2227 

G.1797 

0.1172 

0.0767 

0.0331 

0.0145 

1 

0.30 

1.8148 

1.4885 

1.2223 

1.0047 

0.8268 

0.5618 

0.3835 

0.1811 

0.0871 

1 

0.50 

2.2159 

1.8324 

1.5171 

1.2577 

1.0439 

0.7218 

0.5015 

0.2457 

0.1227 

0.75 

2.3210 

1.9308 

1.6083 

1.3413 

1.1202 

0.7842 

0.5519 

0.2774 

0.1421 

t. 

1.00 

2.3011 

1.9217 

1.6070 

1.3456 

1.1282 

0.7962 

0.5649 

0.2886 

0.1502 

0.0796 

Ld 

1.50 

2.1798 

1.8295 

1.5376 

1.2940 

1.0905 

0.7775 

0.5572 

0.2905 

0.1543 

MEMMi 

2.00 

2.0436 

1.7205 

1.4505 

1.2245 

1.0351 

0.7426 

0.5355 

0.2826 

0.1519 

0.0830 

3.00 

1.7928 

1.5147 

1.2816 

1.0858 

0.9211 

0.6655 

0.4833 

0.2585 

0.1408 

■ 

4.00 

1.5900 

1.3461 

1.1412 

0.9688 

0.8236 

0.5974 

0.2348 

0.1288 

■  fgijS 

5.00 

1.4365 

1.2180 

1.0341 

0.8791 

0.7484 

0.5444 

0.3980 

0.2157 

0.1189 

1 

r  =  5.0 


<t>d 

0.00 

0.25 

0.50 

0.75 

1.00 

1.50 

2.00 

3.00 

4.00 

5.00 

0.10 

0.4777 

0.3835 

0.3080 

0.2476 

0.1991 

0.1289 

0.0837 

0.0066 

0.30 

2.0449 

1.6705 

1.3660 

1.1179 

0.9158 

0.6163 

0.4162 

0.1920 

0.0427 

0.50 

2.5014 

2.0595 

1.6974 

1.4003 

1.1565 

0.7912 

0.5435 

0.2596 

0.1259 

0.0620 

0.75 

2.6237 

2.1725 

1.8008 

1.4944 

1.2414 

0.8595 

0.5975 

0.2924 

0.1454 

0.0734 

td 

1.00 

2.6034 

2.1638 

1.8004 

1.4997 

1.3506 

0.8735 

0.6114 

0.3039 

0.1535 

0.0786 

1.50 

2.4690 

2.0620 

1.7240 

1.4431 

1.2093 

0.8520 

0.6029 

0.3056 

0.1573 

0.0821 

2.00 

2.3163 

1.9403 

1.6272 

1.3661 

1.1482 

0.8139 

0.5793 

0.2971 

0.1547 

0.0817 

3.00 

2.0366 

1.7094 

1.4385 

1.2119 

1.0221 

0.7295 

0.5228 

0.2717 

0.1433 

0.0766 

4.00 

1.8043 

1.5197 

1.2814 

1.0817 

0.9141 

0.6548 

0.4711 

0.2466 

0.1310 

0.0705 

5.00 

1.6307 

1.3754 

1.1613 

0.9817 

0.8307 

0.5968 

0.4304 

0.2265 

0.1209 

0.0654 
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<i>d 

Figure  3.  IIo(t(i,$d)  as  a  function  of  <f>d  with  td  as  parameter  for  r 


0.0  0.5  1.0  1.5  Z.o  2.5  3.0  3.5  1.0  « . 


td 

Figure  4.  flo {tdi4>d)  as  a  function  of  td  with  4>d  as  parameter  for  r 


Appendix  A:  Calculation  of  Q{k,t) 
The  basic  building  block  for  the  calculation  of 


Q(k,t)  =  L~1  f 


—  1  f  \/  S  +  Qiiv/5  "h 


)  =L~1(Q(k,s)) 


is  the  formula 


LT1  ( R{s ))  =  L-1 


1  (*  +  «-/?)» 


r  =m 


s  (s  +  a  +  ^)Jj 

~atI0{pt)  +  (a  —  (3)  f  e~auI0{(3u)  du 
Jo 


H(t ) 


(A.l) 


(A.2) 


In  the  foregoing  expression,  I0  is  the  modified  Bessel  function  of  zero  order,  and  H{t )  is 
the  step  function.  If  we  make  the  identification 


a  — 

ax  +  a2 

-  J_j 

-  +  ± 
\0\ 

2 

in  ’ 

-  02 

=  —  1 

2 

2  \i  1 

we  can  write  (suppressing  the  k  dependence) 

Q(s)  =  (s  +  ai)£(s) 

We  then  have 


(A.3) 


(A  A) 


L-1Q{s)=L~1{(s  +  al)R(s))=Q(t)  , 

Q{t)  =  aiR{t)  + -^R{t)  . 

Substituting  equation  (A.2)  into  equation  (A. 5)  yields 


Q(t)  =  e-«t[pi1(0t)+alo(pt)} 


+  {a2-(32)f  e~auIo{0u )  du  +  6(t) 

Jo 


(A.5) 


(A.6) 
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Appendix  B:  Application  of  Faltung  theorem 


In  the  main  body  of  the  report  we  are  concerned  with  calculations  of  the  form 

1  f°° 

fu{x)  =  —  /  fi{x  -  y)h{y)  dy  (B.l) 

2?r  j_oq 

where  f2{y)  is  given  by 

My)  =  ^2  s{y)  •  (B.2) 

For  computational  purposes  it  may  be  desirable  to  organize  the  computation  of  f\i{x) 
so  that  only  the  function  g(y)  appears  in  the  integration,  and  none  of  its  derivatives.  This 
is  accomplished  in  the  following  way:  We  let 

Ay)  =  |  (b.3) 

so  that  we  can  write 


Mv)  =  ^  g'(y)  ■ 

Substituting  equation  (B.4)  in  equation  (B.l)  gives 

h2  =  hS.Jx(x~y)%  iy 

We  now  integrate  equation  (B.5)  by  parts  to  give 


(B.4) 


(B.5) 


hi  =  jV(y)/i(* -  y)l“„  -  J  g'(y)  — ^  iy 
= "2=  Lj{y) 


df i(x  -  y) 
dy 


(B.6) 


Letting 


we  have 


z  =  x  -  y  , 


(B.7) 
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There  results 


d  =  dfj  dz  _  dfx 
dy  dz  dy  dz 


By  executing  the  same  procedure  once  more  we  obtain  the  desired  result 


Mz)  =  hj-J(y){ 


d2fi{z) 

d2z 


dy 


z=x  —  y 


(B.8) 


(B.9) 


(B.10) 
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